# ---- Caption and Heading of File ---------------------------------------------
#                                                                              #
#        Replication script for creating Figure 6.                             #
#        Not called directly, is a subroutine.                                 #
#                                                                              #
################################################################################
  
  # set-up for Figure 6
  xvalues = seq(0, 1, by=.01)
  empirical_results <- data.frame(pg=xvalues, PG=NA, SE=NA, MOE=NA, lowerbound=NA, upperbound=NA)
  for(i in 0:100)
  {
    # i am showing the uncertainty from the empirical relationship only,
    # not the uncertainty from small sample_n, so sample_n set to very large value
    empirical_results[i+1, 2:6] <- sim.as.jury.stats(sample_pg = xvalues[i+1], sample_n = 1000000, nDraws = 1000)
  }
  
  par(mar=c(3,3,1,.5), mfrow=c(1,2), omi=base::c(0.1,0.1,0.1,0.1), family="serif")
  arrow.length = .05
  horiz_stop = .1
  vert_start = .035
  sample_pg = .64
  sample_n = 354
  MOE_of_pg <- sate:::se.prop(p = sample_pg, n= sample_n) * qnorm(p = .975); MOE_of_pg
  MOE_of_pg <- .05 # asserted value
  
  basic.plot.grid(main="")
  mtext(text="(a) Estimate with Empirical Model", side=3, line=.5, cex=1, font=2)
  # predictions$fit come from aggregate votes and plot results
  lowerbound <- empirical_results$lowerbound-.01
  upperbound <- empirical_results$upperbound+.01
  polygon(x=c(xvalues, rev(xvalues)), y=c(upperbound, rev(lowerbound)), col="gray80", border=NA)
  lines(x=empirical_results$pg, y=empirical_results$PG, col="black", lwd=1, lty=1)
  # this shows CI of P(g) 
  par(xpd=T)
  pg_placement = 0
  points(x=sample_pg, y=pg_placement, pch=16, col="black", cex=.9)
  segments(x0= sample_pg - MOE_of_pg, x1 = sample_pg + MOE_of_pg, 
           y0 = pg_placement, y1 = pg_placement, lty=3, col = "black")
  par(xpd=F)
  # point estimate, not using tace_arrows function because its endpoint based on formal model
  segments(x0 = sample_pg, x1 = sample_pg,
           y0 = vert_start, y1 = sim.as.jury.stats(sample_pg, sample_n)$PG, col="black")
  arrows(y0 = sim.as.jury.stats(sample_pg, sample_n)$PG, y1 = sim.as.jury.stats(sample_pg, sample_n)$PG,
         x0 = sample_pg, x1 = horiz_stop, length = arrow.length, col="black")
  # lower and upper bounds show added uncertainty
  # lower bound
  segments(x0 = sample_pg - MOE_of_pg, x1 = sample_pg - MOE_of_pg,  
           y0 = vert_start, y1 = sim.as.jury.stats(sample_pg - MOE_of_pg, 10000000)$'Lower 5%',  
           col="black", lty=3)
  arrows(y0 = sim.as.jury.stats(sample_pg - MOE_of_pg, 10000000)$'Lower 5%', 
         y1 = sim.as.jury.stats(sample_pg - MOE_of_pg, 10000000)$'Lower 5%', 
         x0 = sample_pg - MOE_of_pg, x1 = horiz_stop,   
         length = arrow.length, col="black", lty=3)
  # upper bound
  segments(x0 = sample_pg + MOE_of_pg, x1 = sample_pg + MOE_of_pg,  
           y0 = vert_start, y1 = sim.as.jury.stats(sample_pg + MOE_of_pg, 10000000)$'Upper 95%',  
           col="black", lty=3)
  arrows(y0 = sim.as.jury.stats(sample_pg + MOE_of_pg, 10000000)$'Upper 95%', 
         y1 = sim.as.jury.stats(sample_pg + MOE_of_pg, 10000000)$'Upper 95%', 
         x0 = sample_pg + MOE_of_pg, x1 = horiz_stop,   
         length = arrow.length, col="black", lty=3)
  # this shows CI of P(G) to left of y-axis
  par(xpd=T)
  horizontal_placement <- -.025
  points(y=sim.as.jury.stats(sample_pg, sample_n)$PG,
         x=horizontal_placement, pch=16, col="black", cex=.9) #
  segments(y0=sim.as.jury.stats(sample_pg - MOE_of_pg, 10000000)$'Lower 5%',  x0 = horizontal_placement,
           y1=sim.as.jury.stats(sample_pg + MOE_of_pg, 10000000)$'Upper 95%', x1 = horizontal_placement, lty=3, col="black")
  par(xpd=F)
  text(x=.01, y=.62, pos=1, label="+/- .115 Margin of Error", cex=.8, srt=90) 
  
  
  basic.plot.grid(main="")
  graphics::mtext(text="(b) Estimate with Formal Model", side=3, line=.5, cex=1, font=2)
  pg_values_scale <- base::seq(0,1,by=.01)
  pG_values_scale <- sate::as.jury.point(pg_values_scale, jury_n=12, pstrikes=0, dstrikes=0, accuracy=.15)
  graphics::lines(x=pg_values_scale, y=pG_values_scale, col="black", lwd=1, lty=1)
  # this shows CI of P(g) below x-axis
  par(xpd=T)
  pg_placement = 0
  points(x=sample_pg, y=pg_placement, pch=16, col="black", cex=.9)
  segments(x0= sample_pg - MOE_of_pg, x1 = sample_pg + MOE_of_pg, 
           y0 = pg_placement, y1 = pg_placement, lty=3, col = "black")
  par(xpd=F)
  trace_arrows(sample_pg = sample_pg, col = "black")
  trace_arrows(sample_pg = sample_pg - MOE_of_pg, col = "black", lty=3)
  trace_arrows(sample_pg = sample_pg + MOE_of_pg, col = "black", lty=3)
  # this shows CI of P(G) to left of y-axis
  par(xpd=T)
  horizontal_placement <- -.025
  points(y=as.jury.point(sample_pg), x=horizontal_placement, pch=16, col="black", cex=.9) #
  segments(y0=as.jury.point(sample_pg - MOE_of_pg),  x0 = horizontal_placement,
           y1=as.jury.point(sample_pg + MOE_of_pg), x1 = horizontal_placement, lty=3, col="black")
  par(xpd=F)
  text(x=.01, y=.62, pos=1, label="+/- .098 Margin of Error", cex=.8, srt=90)
  


